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ABSTRACT 



We suggest a simple fitting formula to represent Comptonized X- and 7-ray spectra 
from a hot (kT e ^ 10 keV), Thomson thick (vr J> 5) hybrid thermal/nonthermal 
plasma in spherical geometry with homogeneous soft-photon injection throughout the 
Comptonizing region. We have used this formula to fit a large data base of Monte-Carlo 
generated photon spectra, and provide correlations between the physical parameters 
of the plasma and the fit parameters of our analytic fit function. These allow us 
to construct Thomson thick Comptonization spectra without performing computer 
intensive Monte Carlo simulations of the high-Tp hybrid-plasma Comptonization 
problem. Our formulae can easily be used in data analysis packages such as XSPEC, 
thus rendering rapid x 2 fitting of such spectra to real data feasible. 



1. Introduction 



Compton upscattering of soft (optical, UV or soft X-ray) radiation by hot (kT e ^ 10 keV) 
plasma is believed to play an important role in the formation of high-energy (hard X-ray and 
7-ray) spectra of many astrophysical objects. Some examples are Galactic X-ray binaries (for 
recent reviews see |Tanaka fc Lewin 1995j |Liang 1998|) , supernova remnants (e.g., [Hillas et al.| 
1998 , Aharonian 1999j ), Seyfert galaxies (for recent reviews see Mushotzky et al. 1993| , |Zdziarskl 
1999), blazars ( |Maraschi et al. 1992 , Dermer et al. 19*92] ), the intergalactic medium ( punyaev &j 
Zeldovich 1980| , Rcphacli 1995), and possibly also 7-ray bursts ( [Liang 1997 , Liang et al. 1999| ). 
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Comptonization in purely thermal plasmas is fairly well understood, and analytical solutions 
for the emerging hard X-ray and 7-ray spectra for a wide range of Thomson depths and plasma 
temperatures have been developed ( [Sunyaev Titarchuk 1980| , [Titarchuk 1994 , Titarchuk & 



Ljubarskij 1995] , [Hua fc Titarchuk 1995| ), However, the currently known solutions can not be 



applied to hybrid thermal/nonthermal or even purely nonthermal plasmas, which are believed to 
be the primary source of hard X-ray radiation in many astrophysical objects. For example, Li, 
Kusunose, & Liang (1996) have shown that the high-energy spectrum of Cyg X-l can be plausibly 
reproduced by Comptonization in a quasi-thermal plasma with a small fraction of electrons in a 
nonthermal tail. Thermal/nonthermal hybrid models have also been discussed in the context of the 
high-energy spectra of Seyfert galaxies, e.g. NGC 4151 ( |Zdziarski et al. 1994 , [Zdziarski, Johnson 



fc Magdziarz 1996| ). The development of such nonthermal tails has been confirmed by Bottcher 



& Liang (2001) by detailed Monte-Car lo/Fokker-Planck simulations of radiation transport and 
electron dynamics in magnetized hot plasmas, e.g. in accretion flows onto compact objects. 

We need to point out that the combination of soft- to medium-energy X-ray spectra (typically 
integrated over a few ksec) with hard X-ray / soft 7-ray OSSE spectra (typically integrated 
over hundreds of ksec or even several days) could introduce artificial spectral signatures due to 
the superposition of thermal Comptonization spectra from different emission regions, possibly 
at different times. However, the detailed Fokker-Planck simulations of the dynamics of electron 
acceleration and cooling in coronal plasmas in the vicinity of accreting compact objects (e.g., |L1| 
et al. 1996 , Bottcher fc Liang 200 1| ) clearly demonstrate that nonthermal tails in the electron 



distributions are very likely to be produced in such environments. It is thus desirable to have 
simple, analytical expressions for the Comptonization spectra arising in this case. 

For the general case of arbitrary Thomson depth, average electron energy (i.e. temperature), 
nonthermal fraction of electrons in the hybrid plasma, and nonthermal electron spectral 
parameters, the most reliable and often fastest method for computing Comptonization spectra 
is the Monte-Carlo method ( Pozdnyakov et al. 1983j ). However, Monte-Carlo simulations of 



Comptonization problems become extremely time consuming in the case of very high Thomson 
depth, tt, since the number of scatterings that need to be simulated, increases oc r\. To date, 
there are no simple, analytic approximations to the problem of optically thick Comptonization in 
hot, thermal/nonthermal hybrid plasmas. This has made x 2 minimization fitting of such models 
to real data infeasible until now. 

Using our linear Monte-Carlo Comptonization codes, we have built up a data base of over 
300 Comptonization spectra from hot, Thomson thick, thermal/nonthermal hybrid plasmas 
with different Thomson depths, temperatures, nonthermal electron fractions, spectral indices 
and maximum electron energies of nonthermal electrons. In this paper, we propose a simple 
analytic representation, consisting of an exponentially cut-off power-law plus a smoothly connected 
double-power-law ( |Band et al. 1993| ) to fit all individual spectra in our data base. We find 
correlations between the physical parameters of the Monte-Carlo calculations and the parameters 
of our analytical representation, which can be used to construct Comptonization spectra of 
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high-rx plasmas with arbitrary nonthermal electron fractions, without having to perform 
time-consuming Monte-Carlo simulations. Our results are specifically optimized for very high 
Thomson depths, ^ 5, and are valid for electron temperatures kT e ;> 10 keV, since we assume 
the Comptonizing plasma to be fully ionized. For lower temperatures, effects of photoelectric 
absorption, recombination, and other atomic processes, which have been neglected in the current 
work, must be included. 

2. Model setup and analytic representation of the Comptonization spectra 

The physical scenario underlying the Comptonization problem at hand, consists of a spherical, 
homogeneous region of hot (kT e ;> 10 keV), thermal/nonthermal hybrid plasma with radial 
Thomson depth tt- The thermal/nonthermal hybrid distribution function of electron energies is 
given by 



Here, the normalization factors N t h and N n t and the transition energy 71 are determined by the 
normalization, / e (7) d-f = 1, by the requirement of continuity of the distribution function at 
71, and by the fraction a mxw i of thermal particles, i. e. f^ 1 f e {l) d-7 = flmxwi- 

Soft photons with characteristic energy £ s < 1 keV are injected uniformly throughout the 
Comptonizing region (we assume that they have a thermal blackbody spectrum with kT r <C 1 keV; 
the specific shape of the soft photon distribution is irrelevant for multiple-Compton-scattering 
problems). The Comptonizing plasma is assumed to be fully ionized. The Comptonization spectra 
are calculated with our Monte-Carlo codes, which use the full Klein-Nishina cross section for 
Compton scattering. 77 absorption, pair production, and other induced processes (e.g., free- free 
absorption and induced Compton scattering events) have been neglected in our Monte-Carlo 
simulations. 

We have created a data base of over 300 simulations, spanning values of tt < 28, 
5keV < kT e < 200 keV, 0.2 < a mxwl < 1, 2.5 < p < 6.5, and 10 < 72 < 10 3 . 

We here propose a fit function to optically thick Comptonization spectra from 
thermal/nonthermal hybrid plasmas as the sum of an exponentially cut-off power-law plus 
a Band GRB function (defined here, for convenience, in energy flux units to match the units of 
the output spectra of our Monte-Carlo simulations): 




(1) 



F E (E) = N PL E- r e- E/Eo 



- 4 - 



f E 5 e- E / Eo for E <{8- e)E 
+ A- I (2) 
[ B E e for E > (5 - e)E 

with B = {\8 — e\Eo) S ~ <L e e ~ s and E in keV. .Fg is the energy flux in units of ergs cm -2 s _1 keV -1 . 

Apart from the overall normalization, there are 5 parameters in this fit function: the 
spectral indices T, 5, and e, the turnover energy Eq (corresponding to a vF v peak energy of 
-E'pk = [1 + 5]Eq), and the relative normalization of the two components, A/Npi,. We have 
developed a code using a combination of a coarse grid in the 6-dimensional parameter space with 
a x 2 forward folding method to fit this function to our simulated Compton spectra. Fig. |l] shows 
two typical examples of the resulting fits. 



3. Correlations between physical and fit parameters 

We have fitted the function (^) to the energy range 1 keV < E < 10 MeV of all simulated 
Comptonization spectra in our data base. In order to be able to use our fit function for physically 
meaningful fitting, one has to establish a unique correlation between the physical parameters of 
the Monte-Carlo simulation and the best-fit parameters of our fit function. First of all, we note 
that in the low-energy range, dominated by the power-law part of the spectrum, Fe oc E~ r , 
the influence of the nonthermal population is negligible, and we can use the standard result for 
thermal plasmas in spherical geometry in the case of homogeneous photon injection throughout 
the Comptonizing region: 




4vr 2 



27 G e (r T + 2/3) 2 J 

( Sunyaev Sz Titarchuk 198C| ) where G e = kT e /(m e c 2 ). This relation is valid in the limit of 



(3) 



high Thomson depth, tt ^ 5, and moderate electron temperature, kT e <j 250 keV ( |Sunyaev~& 



Titarchuk 1980, Pozdnyakov et al. 1983j ). A fully analytical solution for the power-law index T 



for all values of tt and e has been derived by Titarchuk & Ljubarskij (1995). In the parameter 
range in which our analytical representation is applicable, the simple expression (|3|) is sufficiently 
accurate. 

A second, obvious correlation exists between the vF v peak of an optically thick Comptonization 
spectrum and the average particle energy. In the thermal case, the Wien spectrum peaks at 
E-pk « 3 kT e . However, near the vF v peak of the spectrum, the influence of the nonthermal particle 
population can become important. Generalizing the E p ^ vs. kT e correlation, one would expect 
that Epk ~ 2(7 — l)m e c 2 , where the brackets denote the ensemble average. We find that this 
generally provides a reasonable fit to the simulated E p ^ values. However, in some cases (especially 
for hard power-law tails, p <^ 3), the onset of the nonthermal population (at 71) is at too high an 
energy and the normalization of the nonthermal population is too small for nonthermal electrons 
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to make a strong contribution to photons scattered into the Wien peak. We find significant 
deviations from the above E p k vs. (7 — 1) correlation for values of 71 > 6 e + 1. An approximate 
correction for this effect can be found in the following empirical correlation which provides a 
reasonable description of the simulated E p ^ values: 

( 1.022 (7 - 1} MeV for 71 < 6 G e + 1 

E pk = \ _ (4) 

°mxwl(7-l>th + (l-a mX wl)(7-l)nt2-//8 M y f g @ 1 

I. a mxw i+(l-o mxw i)2--r/ 8 n - 

where / = 2(71 — l)/@ e . The terms (7 — l) t h and (7 — l) n t denote the average kinetic energies 
(normalized to m e c 2 ) of electrons in the thermal and nonthermal portion of the electron 
distribution, respectively. The E p k vs. 7 — 1 correlation for those cases with 71 < 6 O e + 1 is 
illustrated in Fig. We note that in most physically relevant cases, we expect the power-law 
index of nonthermal electrons to be p ^ 3 since we do not expect a positive vF v spectral slope at 
high energies due to the nonthermal portion of the hybrid electron spectra. In the case p ^ 3, the 
upper branch of Eq. (||) always yields satisfactory results. 

The value of (7 — 1) depends on the temperature of the thermal population and the parameters 
of the nonthermal population in a non-trivial way. We find that (7 — 1) can be parametrized 
rather accurately by the fitting function 



(7-I) kK x ®^ +K 2 e e ^ K2 . (5) 

To facilitate the use of our fitting formulae, we list in Table [l] the relevant fit parameters K\, 
K2, Ki, and K2 for a representative value of 72 = 10 3 . The dependence of (7 — 1) on the thermal 
plasma temperature and the nonthermal power-law index p is illustrated in Fig. ||. 

The connections of the remaining fit parameters, e, 5, and A/N-pi, with the physical parameters 
of the Comptonizing plasma, are less obvious, and the correlations which we have found on the 
basis of the fit results to our data base, are purely empirical. In an optically thin (tt <C 1), purely 
nonthermal plasma, the high-energy spectral index e would be related to the non-thermal particle 
index p through e = — (p — l)/2 (e.g., Rybicki fc Lightman 1979| ), However, the presence of the 



thermal population as well as the effect of multiple Compton scatterings lead to a systematic 
distortion of this spectral shape. We find a correction to e due to these effects, which can be 
parametrized as 



with 



- - (a e + b e e e ) (6) 



a e = 0.78 e- (L49a — 1)3 (7) 

-0.7 
mxwl 



be = 3.25 a- ^ (8) 
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Two examples of the correlation (||) are shown in Fig. |||. 

The low-energy power-law slope of the Wien hump, added to the low-energy power-law (the 
-Fe oc E s part of our fitting function), approaches 5 — > 2 for very high Thomson depth. We find 
that the correlation between 5 and tt can be well parametrized by a functional form 

(9) 



[ro(0e )],(ee) +r ^) 

which is virtually independent of a mxw i and the specifications of the nonthermal population. The 
parameters in Eq. |9] depend on the thermal electron temperature as 



r (G e ) = 5.6 9°' 13 (10) 
r?(9 e ) = 2 (ln[200e e ]) L7 (11) 

Finally, we expect the ratio of power in the Wien hump (Band function section of our fit 
function) to the power in the low-energy power-law part of the Comptonization spectrum to 
be positively correlated with the plasma temperature G e and Thomson depth tt. The above 
mentioned power ratio can be expressed as the term 

C A/NpL ^(A/N PL )E s p fe-( 5+1 l (12) 

To first order, this power ratio should not strongly depend on the specifications of the non-thermal 
component. Fig. || shows some examples of the correlation between Ca/n pl an d tt for various 
plasma temperatures, where we have neglected any dependence on a mxw i and p. We find that this 
correlation can be parametrized as 

C'a/n pl = C {^) C + Ci (13) 

with 



Co = 16 G" 1 (14) 
Ci = 5 (15) 
a c = e°- 05 (16) 

However, comparing the analytic spectra calculated using Eq. (^) with the Monte-Carlo generated 
spectra, we do find significant a mxw i-dependent deviations from the best-fit values of A/Npi,. We 
thus introduce a correction C%/jv pl such that 
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Ca/n pl — C'a/n pl • C" A /N PL (17) 

and find this correction term as 

C\ /JVpL = (2.1 • 10- 3 9- 2 - 3 + 1.8) r- A(0e) (18) 

with 

A(e e ) = 2.8-io- 2 e-°- 95 . (19) 



Since we know E p k, 5 and T from Eqs. (Q), @ and (^), we can easily invert Eq. (12) to find the 
fit parameter A/Np^. 



4. Summary and Conclusions 

We have developed an analytical representation of Comptonization spectra from optically 
thick, hot, thermal/nonthermal hybrid plasmas in spherical geometry with homogeneous 
soft-photon injection. Starting from the physical parameters of the problem at hand, one can 
easily use Eqs. (H|), (||), @, (f|), and ( |l3|) to construct the resulting photon spectrum. Thus, our 
representation can be used in x 2 minimization software packages to scan through the physical 
parameter space of the thermal/nonthermal hybrid plasma and fit the emerging Comptonization 
spectra to observed photon spectra. Fig. |6] shows two representative examples of the analytical 
spectra, using the correlations presented in the previous section, and the output spectra from the 
corresponding Monte-Carlo simulations. 

Our results are valid for electron temperatures lOkeV <J kT e ^ 150 keV, Thomson depths 
tt ^ 5, Maxwellian electron fractions 0.5 <; a mxw i <^ 0.95 (for purely thermal plasmas, a mxw i = 1, 
the analytical solutions of Hua & Titarchuk [1995] may be used), nonthermal electron spectral 
indices p ^ 3, and seed photon energies £ s « 1 keV. The results have been obtained for 
homogeneous photon injection in spherical geometry and should not be used for other geometries. 
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E [keV] 

Fig. 1. — Two examples of fits of our model function (solid and dashed curves) to the Monte- 
Carlo simulated energy spectrum (filled symbols). Parameters for the Monte-Carlo simulation 
(filled circles) are: r T = 16, kT e = 80 keV, a mxw i = 0.5, p = 3.5, 72 = 10 3 , kT s = 0.5 eV. 
Parameters of the analytic fit (dashed curve) are: T = —0.035, S = 1.98, e = —2.27, Eq = 196 keV, 
A/N p i = 8.84 x 10~ 4 . Parameters for the second MC simulation (filled squares) are: tt = 28, 
kT e = 40 keV, a mxw i = 0.7, p = 3.5, 72 = 10 3 , kT s = 0.5 eV. Parameters of the analytic fit (solid 
curve): T = -0.029, 5 = 2.01, e = -1.90, E = 149 keV, A/N pl = 1.87 x 10~ 3 . 
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r 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 




<y- 1> 

Fig. 2. — Correlation between the average electron energy in the plasma, (7 — 1), and the vF v peak 
energy, E p ^, compared to the expected correlation, = 1022(7 — 1) keV (solid line). The best-fit 
linear slope is 973 keV (dot-dashed line). 
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Fig. 3. — Dependence of (7 — 1) on plasma temperature and nonthermal spectral index p; example 
for a mxw i = 0.6, 72 = 10 3 . From top to bottom, the nonthermal spectral indices are p = 2.1, 2.23, 
2.35, 2.48, 2.61, 2.73, 2.86, 3.11, 3.37, 3.62, and 3.87. 




Fig. 4. — Correlation between the correction to the high-energy spectral index e and the plasma 
temperature for two different values of the Maxwellian plasma fraction a mxw i. 
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Fig. 5. — Empirical correlation between the relative power in the power-law to the Wien hump 
(Band function part of the photon spectrum) and the Thomson depth tt for various electron 
temperatures (with various values of a mxw l each, which does not appear to have a significant 
influence). The lines are the best-fit power-laws. 
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E [keV] 

Fig. 6. — Two typical examples of the analytical representation for Thomson thick Comptonization 
spectra, using Eqs. (0), (0), (0), (H), and ( |l3| ) (solid and dot-dashed curves), compared to the 
Monte-Carlo simulated output spectra (symbols; same simulations as shownin Fig. |l|). In both 
simulations, p = 3.5 and 72 = 10 3 was used. 
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Table 1: Fit parameters to compute (7 — 1) as a function of the plasma parameters using Eq. [5|. 
72 = 10 3 has been used. 
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